function [A, dAdpar] = T1FAVMagnitude( par, echoes)

% computes the magnitude M with given parameters and inversion times
% echoes = [n x 1] inversion times
% par = the 3 element parametervector [A, B, R1]'
%
% The datamodel is given by
% SI(i) = M(1-E)*sin(echoes(i)) / (1-E*cos(echoes(i)))
%
% Created by Henk Smit, EMC, 01-2011 based on the work by Dirk Poot, University of Antwerp, 13-8-2007

if size(par,1) ~= 2 %HENK
    error('incorrect parameter vector');
end;
numtr = size(par,2);
numgr = size(echoes,1);

SI = zeros(numgr,numtr);
dSIdM=SI;
dSIdE=SI;

for k=1:numtr
    M=par(1,k);
    E=par(2,k);
    angles=echoes(:,1);
    dSIdM(:,k) = ((1-E)*sin(angles))./(1-E*cos(angles));
    dSIdE(:,k) = ((sin(angles).*cos(angles).*M.*(1-E)) - (M.*sin(angles).*(1-E.*cos(angles))))./ ( (1 - (E.*cos(angles))).^2);
    SI = M*dSIdM;
    A=SI;
end;

if nargout>1
    dAdpar=[dSIdM dSIdE]; %HENK 
end;
